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ABSTRACT 

We present accurate time delays for the quadruply imaged quasar HE 0435-1223. The delays were measured from 575 independent 
photometric points obtained in the R-band between January 2004 and March 2010. With seven years of data, we clearly show that 
quasar image A is affected by strong microlensing variations and that the time delays are best expressed relative to quasar image B. 
We measured At BC = 7.8 + 0.8 days, At BD = -6.5 + 0.7 days and Af CD = -14.3 + 0.8 days. We spacially deconvolved HST NICMOS2 
F160W images to derive accurate astrometry of the quasar images and to infer the light profile of the lensing galaxy. We combined 
these images with a stellar population fitting of a deep VLT spectrum of the lensing galaxy to estimate the baryonic fraction, f b , 
in the Einstein radius. We measured f b = 0.65^-Jq if the lensing galaxy has a Salpeter IMF and f b = 0.45 ^-j^ if it has a Kroupa 
IMF. The spectrum also allowed us to estimate the velocity dispersion of the lensing galaxy, cr ap = 222 ± 34 km s _1 . We used f b and 
cr ap to constrain an analytical model of the lensing galaxy composed of an Hernquist plus generalized NFW profile. We solved the 
Jeans equations numerically for the model and explored the parameter space under the additional requirement that the model must 
predict the correct astrometry for the quasar images. Given the current error bars on f b and cr ap , we did not constrain H yet with high 
accuracy, i.e., we found a broad range of models with^ 2 < 1. However, narrowing this range is possible, provided a better velocity 
dispersion measurement becomes available. In addition, increasing the depth of the current HST imaging data of HE 0435-1223 will 
allow us to combine our constraints with lens reconstruction techniques that make use of the full Einstein ring that is visible in this 
object. 

Key words. Gravitational lensing: time delay - Cosmology: Hubble constant - Galaxies: quasar: individual (HE 0435-1223) 



1. Introduction 

To determine the expansion rate of the Universe, Hq, accurately, 
it is important to scale the extragalactic distance ladder and to 



* Based on observations made with the 1.2m Euler Swiss Telescope, 
the 1.5m telescope of Maidanak Observatory in Uzbekistan, and with 
the 1.2m Mercator Telescope, operated on the island of La Palma by 
the Flemish Community, at the Spanish Observatorio del Roque de los 
Muchachos of the Instituto de Astroffsica de Canarias. The NASA/ESA 
Hubble Space Telescope data was obtained from the data archive at 
the Space Telescope Science Institute, which is operated by AURA, 
the Association of Universities for Research in Astronomy, Inc., under 
NASA contract NAS-5-26555. 

** Research Fellow, Belgian National Fund for Scientific Research 
(FNRS) 

*** Alexander von Humboldt fellow 



measure w(z), the redshift evolution of the dark energy equation 
of state parameter (e.g., Frieman 'et al.|2008] ). 

Some of the most popular methods in use to measure #0 
have recently been reviewed (Freedman & Madore 2010). The 
huge observational and theoretical efforts invested in these mea- 
surements have led to random errors between 3% and 10%, de- 
pending on the methods. However, all these methods rely at 
some level on each other and are in the end all based on the same 
local standard candles. In addition, the accuracy required on H 
to measure w(z) a nd H(z) is of the order of 1% ( |Hu|2005t |Riess| 
etal. 2009, 20TT]). It is therefore of interest (i) to explore meth- 
ods that are fully independent of any standard candle and (ii) to 
combine the different methods to further reduce the current error 
bar on Hq. 
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Fig. 1. Part of the field of view of the 1.2 m Swiss Euler telescope, with HE 0435-1223 visible in the centre. The four PSF stars 
used for deconvolution purposes and the two reference stars used to carry out the flux calibration are indicated. 



Strong gravitational lensing of quasars and the so-called 
"time delay method" in multiply imaged quasars is independent 
of the traditional standard candles ( Refsdal 1964) and is based 
on well understood physics: general relativity. However, it re- 
quires to measure the time delays from long-term photometric 
monitoring of the many lensed quasars, which has long been 
a serious observational limitation. While the attempts to mea- 
sure accurate time delays have been numerous, only few quasars 
have been measured with an accuracy close to the percent (e.g., 
Goicoechea 2002; Fassnacht et al. 2002; Vuissoz et al. 2008 ). In 
addition, the lens model necessary to convert time delays into 
Hq often remains poorly constrained and hampers a breaking 
of the degeneracies between the properties of the lensing galaxy 
and Hq. In recent years, successful attempts have been made to 
constrain the lens model as much as possible to break these de- 
generacies ( |Suyu et al.|2009l[20T0| ). 

COSMOGRAIL, the COSmological MOnitoring of 
GRAvItational Lenses, aims both at measuring precise time 
delays for a large sample of strongly lensed quasars, and at 
obtaining and using all necessary observations to constrain the 
lens models. The present paper describes the COSMOGRAIL 
results for the quadruply imaged quasar HE 0435-1223, using 
deep VLT spectra and deconvolved HST images. 

HE 0435-1223 (a(2000) : 04 ft 38 mm 14.9 sec; £(2000) = 



ing the Hamburg/ESO Survey (HES) for bright quasars in the 
Southern Hemisphere. It w as identified two yea rs later as a 
quadruply imaged quasar by Wisotzki et al. (2002). The redshift 
of the source is z s = 1.689 ( [Wisotzki et al. 12000) an d that of the 
lens is z { = 0.4546 ± 0.0002 (Morga netaL|2005) . The quasar 
shows evidence for intrinsic variability, which makes it a good 
candidate for determining the time delays between the different 
images. The local environment of the lensing galaxy has been 
studied in detail by Morgan et al. (2005 ) using the Hubble Space 
Telescope (HST) Advanced Camera for Surveys (ACS) and by 
Momcheva ( 2009), who has found that the lensing galaxy lies in 
a group of at least 1 1 members. The velocity dispersion of this 
group is cr ~ 496 km s" 1 . 



Analytical lens models of HE 0435-1223 are given by 
Kochanek et al.| (|2006|), who have also measured time delays 



-12°17'14'/4) was discovered by Wisotzki et al. (2000b dur- 



from two years of optical monitoring: At ad = -14.37^ ^ days, 
At AB = -8.00+g-g days, and At AC = -2.10+°-™ days. For a fixed 
Hq = 72 ± 7 km s _1 Mpc _1 they found that the lensing galaxy 
must have a rising rotation curve at the position of the lensed 
images and a non-constant mass-to-light ratio. Moreover, high 
dark matter surface densities are required in the lens halo. New 
monitoring data of Blackburne & Kochanek] ( |20 10) analysed us- 
ing a physically motivated representation of microlensing give 
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Fig. 2. Result from the simultaneous deconvolution of the ground-based frames. G is the lensing galaxy and G22 (Morgan et al. 
2005| ) is its closest neighbour on the plane of the sky. The grey scale in the deconvolved image is set to display all light level above 
3 x cr sky- The FWHM resolution of the deconvolved image is Of! 34. 

Table 1. Summary of the optical monitoring data. The temporal sampling is the mean number of days between two consecutive 
observations. 



Telescope 


Camera 


FoV 


Pixel 


Period of observation 


#obs. 


Exp. time 


Seeing 


Sampling 


Euler 


C2 


li'x ir 


a/344 


Jan 2004 - Mar 2010 


301 


5 x 360s 


1737 


6 days 


Mercator 


MEROPE 


6.5'x 6.5' 


07 190 


Sep 2004 - Dec 2008 


104 


5 x 360s 


1759 


1 1 days 


Maidanak 


SITE 


8.9'x 3.5' 


a/266 


Oct 2004 - Jul 2006 


26 


10 x 180s 


1731 


16 days 


Maidanak 


SI 


is. i'x i8.r 


a/266 


Aug 2006 - Jan 2007 


8 


6 x 300s 


1731 


16 days 


SMARTS 


ANDICAM 


10'x 10' 


a/300 


Aug 2003 - Apr 2005 


136 


3 x 300s 


< 1780 


4 days 


TOTAL 








Aug 2003 - Mar 2010 


575 


242.5 h 




3.2 days 



time delays compatible with those of |Kochanek et al.| (2006), 
although these authors do not provide their measured values. 

2. Photometric monitoring 

2.1. Opticalimaging 

HE 0435-1223 was monitored during more than six years, from 
January 2004 to March 2010, through the R filter, using three 
different telescopes: the Swiss 1.2m Euler telescope located on 
the ESO La Silla site (Chile), the Belgian-Swiss 1.2m Mercator 
telescope located at the Roque de Los Muchachos Observatory, 
La Palma, Canary Island (Spain), and the 1.5m telescope located 
at the Maidanak Observatory (Uzbekistan). In addition we also 
use 136 epo chs from the two-year long monitoring of Kochanek 
|eTaIl ( [20n6l ), from August 2003 to April 2005, obtained with the 
ANDICAM camera mounted on the 1.3m Small and Moderate 
Aperture Research Telescope System (SMARTS) located at the 
Cerro Tololo Inter- American Observatory (CTIO) in Chile. A 
summary of the observations is given in Table [T] Note that we 
used the published photometry for the SMARTS data, i.e., we 
did not reprocess the original data frames with our own photo- 
metric pipeline. 

2.2. Image processing and deconvolution photometry 

The data from the three telescopes used by the COSMOGRAIL 
collaboration are an alysed with the semi - automated reduction 
pipeline described in Vuissoz et al.| ([2007). The main challenge 
was that we had to assemble data from different telescopes: each 
camera has a different size, resolution and orientation on the 
plane of sky. 



The pre-reduction for each observing epoch consists of 
flat-fielding using master sky-flats. The Euler image with the 
best seeing was taken as the reference frame to register all 
other frames. This reference frame was taken on the night of 
November 11, 2005 and has a seeing of 0782. Two reference 
stars were then chosen in the field of view (Fig. [TJ to compute 
the geometrical transformations between the images. This trans- 
formation involves a spatial scaling and a rotation. The reference 
stars were also used to compute the relative photometric scal- 
ing between the frames t aken at different epo chs. Eventually, the 
L.A. Cosmic algorithm ( van Dokkum 2001) was applied sepa- 
rately to every frame to remove cosmic rays. All images were 
checked visually to make sure that no pixel was removed inap- 
propriately, especially in the frames with good seeing. 

The photometric measurements were carried out using "de- 
convolution photometry" with the MCS deconvolution algorithm 
(Maga in~et al.||1998] ). This software has been successfully ap- 
plied to a variety of astrophysic al problems ranging from grav- 
itationally lensed quasars (e.g.,|Burud et al.||2000 



2002 ) to the 
2008 ), or to 



study of quasar host galaxies (e.g., |Letawe et al.] 
the search for extraso lar planets using the transit technique (e.g., 
Gill on et al. 2007a|b| ). Image deconvolution requires accurate 



knowledge of the instrumental and atmospheric point spread 
function (PSF). The latter was computed for each frame from 
the four stars labelled PSF 1-4 in Fig. [I] These stars are from 
to 1 magnitude brighter than the quasar images in HE 0435-1223 
and are located within 2' from the centre of the field, which min- 
imizes PSF distortions. 

Because it does not attempt to achieve an infinitely high spa- 
tial resolution, the MCS algorithm produces deconvolved im- 
ages that are always compatible with the sampling theorem. This 
avoids deconvolution artefacts and allowed us to carry out ac- 
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Fig. 3. R-band light curves of the four lensed images of HE 0435-1223 from December 2003 to April 2010. The magnitudes are 
given in relative units as a function of the Heliocentric Julian Day (HJD), along with their total l<x error bars. 
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Fig. 4. Left: Combination of the four original HST/NIC2 F160W frames of HE 0435-1223. The field of view is 9 x 9 arcsec. Middle: 
Deconvolved image, where the lensing galaxy is modelled as a de Vaucouleurs profile (see text). The nearest galaxy on the plane of 
the sky, G22, is also indicated. Right: Residual map in units of the noise. The colour scale ranges from -4<x (white) to +4<x (black). 



Table 2. Relative astrometry of HE 0435-1223as derived from the simultaneous deconvolution of all NIC2 frames. The l<x error bars 
are the internal errors after deconvolution. Additional 2-mas systematic errors must be added to these (see text). The magnitudes 
are in the V ega system. For comparison we show the results from Morgan et al. ( 2005 ) using HST/ACS images and from Kochanek 
etaL| p006| ) using HST/NIC2 images. 



This work 



Morgan et j 



Kochanek et al. (2006 



ID 


Aa (") 


AS (") 


Mag (F160W) 


Aa(") 


AS (") 


Aor (") 


A* n 


A 


0. 


0. 


17.20 ± 0.01 


0. 


0. 


0. 


0. 


B 


-1.4743 + 0.0004 


+0.5518 + 0.0006 


17.69 ± 0.01 


-1.477 ±0.002 


+0.553 ± 0.002 


-1.476 ±0.003 


+0.553 ± 0.001 


C 


-2.4664 ± 0.0003 


-0.6022 ± 0.0013 


17.69 ± 0.02 


-2.469 ± 0.002 


-0.603 ± 0.002 


-2.467 ± 0.002 


-0.603 ± 0.004 


D 


-0.9378 ± 0.0005 


-1.6160 ± 0.0006 


17.95 ± 0.01 


-0.938 ± 0.002 


-1.615 ±0.002 


-0.939 ± 0.002 


-1.614 ±0.001 


G 


-1.1706 ±0.0030 


-0.5665 ± 0.0004 


16.20 + 0.12 


-1.169 ±0.002 


-0.572 ± 0.002 


-1.165 ±0.002 


-0.573 ± 0.002 



curate photometry over the entire field of view. Moreover, the 
deconvolved image was computed as the sum of extended nu- 
merical structures and of analytical point sources whose shape is 
chosen to be symmetrical Gaussians. In the case of gravitation- 
ally lensed quasars, the numerical channel of this decomposition 
contains the lensing galaxy. The photometry and astrometry of 
the quasar images were returned as a list of intensities and posi- 
tions of Gaussian deconvolved profiles. Finally, the deconvolved 
image can be computed on a grid of pixels of arbitrary size. In 
the present work, the pixel size in the deconvolved frames is half 
the pixel size of the Euler data, i.e., Off 112. The spatial resolu- 
tion in the deconvolved frames is two pixels full-width-at-half- 
maximum (FWHM), i.e., a/35. 

With the MCS software, dithered images of a given target 
can be "simultaneously deconvolved" and combined into a sin- 
gle deep and sharp frame that matches the whole dataset at once, 
given the PSFs and the noise maps of the individual frames. In 
doing this, the intensities of the point sources are allowed to vary 
from one frame to the next while the smooth background, which 
includes the lensing galaxy, is held constant in all frames. The 
result of the process is shown in Fig. [2) where the point sources 
are labelled as in |Wisotzki et a l. (2002). Prior information on the 
object to be deconvolved can be used to achieve the best possible 
results. In the case of HE 0435-1223 the relative positions of the 
point sources are fixed to the HST astrometry obtained in Sect. 

EJ 

Fig. [3] shows the deconvolution light curves obtained for each 
quasar image of HE 0435-1223, where the lcr error bars account 
both for the statistical and systematic errors. The statistical part 
of the error was taken as the dispersion between the photometric 



Table 3. Shape parameters for the lensing galaxy in HE 0435- 
1223. The position angle (PA) is measured positive east of north. 
The lcr error bars (internal errors) are given in parenthesis. 



PA(°) 


Ellipticity 


a e ff (") b eff (") 


r eff (") 


174.8 (1.7) 


0.09 (0.01) 


1.57 (0.09) 1.43 (0.08) 


1.50 (0.08) 



points taken during each night. The systematic errors were esti- 
mated by carrying out the simultaneous deconvolution of refer- 
ence stars in the vicinity of HE 0435-1223. 

Finally, a small scaling factor was applied to the light 
curves of all telescopes , including the published light curves of 
Kochane k et al. ( |2006| ), to match the Euler photometry. These 
shifts are all smaller than 0.03 mag. 



3. HST NICMOS2 imaging 

We used deep near-IR HST images of HE 0435-1223 to de- 
rive the best possible relative astrometry between the quasar 
images and the lensing galaxy and to constrain the light 
distribution in the lensing galaxy. The data are part of the 
CASTLES project (Cf a- Arizona Space Telescope LEns Survey) 
and were acquired in October 2004 (PI: C.S. Kochanek) with 
the camera 2 of NICMOS, the Near-Infrared Camera and Multi- 
Object Spectrometer. They consist of four dithered frames taken 
through the F160W filter (H-band) in the MULTIACCUM mode 
with 19 samples and calibrated by CALNICA, the HST image 
reduction pipeline. The total exposure time amounts to approxi- 
mately 44 minutes and the pixel scale is 07 075652. 
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The MCS deconvolution algorithm was used to combine the 
four NIC2 frames into a deep sharp IR image. We followed 
the iterative technique de scribed in |Chantry et al.| ([2010 ) and 
Chantry & Magain] f2007| ), which allowed us to build a PSF in 
the absence of a stellar image in the field of view. The method 
can be summarised as fo llows. First, we es timated the PSF us- 
ing Tiny Tim software ( |Krist & Hook|2004 ) and carried out the 
simultaneous deconvolution of the four F160W frames using a 
modified version of the MCS software (Magain et al. 2007]). This 
produces a first approximation of the extended channel of the de- 
convolved image, i.e., the lensing galaxy and the lensed quasar 
host galaxy. We reconvolved the latter by the PSF and subtracted 
it from the original data. A new estimate of the PSF was built on 
the new image that now contained only the quasar images. The 
process was repeate d until the residua l image was satisfactory 
(for more details see Chantry et al.|2010| ). Fig. [4] shows the result. 
In this image the pixel size is half that of the original data and 
the resolution 07075 (FWHM), unveiling an almost full Einstein 
ring. 

In the final deconvolved image, the lensing galaxy was mod- 
elled analytically rather than numerically to minimise the num- 
ber of degrees of freedom. We found that the best-fit profile 
is an elliptical de Vaucouleurs with the parameters as given in 
Table [3] The astrometry of the quasar images relative to the 
lensing galaxy, corrected for the known distortions of the NIC2 
camera and for the difference of pixel scale between the x and 
y directions is summarised in Table [2] Based on our previ- 
ous work usin g deconvolution of NICMOS images ( Chantry & 
Magain 2007 ), we estimate that the total error bars, accounting 
for residual correction of the distortions amounts to 2 mas. Our 
resu lts agree well with previous measurements fr om HST/ ACS 
(Mo rgan et aL||2005] ) or HST/NIC2 imaging ( [Kochanek eTaL 
2006), also shown for comparison in Table [2] 



4. Time delay measurement 

4.1. Curve shifting method 

Our method to measure the time delays is based on the disper- 
sion technique of Pe lt et al. ( 1996 ): the light curves are shifted in 
time and in magnitude to minimise a global dispersion function. 
In addition, the light curves are distorted on long time scales 
to account for slow microlensing variations. This was made by 
adding low-order polynomials to either the full curves or to spe- 
cific observing seasons. 



Pelt et al. (T996J) has defined several dispersion statistics be- 



tween pairs of light curves. We i mplemented a dis persion esti- 
mate similar to D\ (see eq. 8 of IPelt et al. 1996), which per- 
formed a linear interpolation between points ol one of the curves 
over a maximum range of 30 days. In the case of four light 
curves, we defined a total dispersion that is the sum of the disper- 
sions computed using the 12 possible permutations of two curves 
among four. Each pair was considered twice so to avoid the ar- 
bitrary choice of a reference light curve. The photometric error 
bars were taken into account to weight the influence of the data 
points in the dispersion. We then minimised the total dispersion 
by modifying the time delays and the microlensing polynomials. 

4.2. Microlensing and influence on the time delay 

Simulated light curves that mimic the observed data were used to 
estimate the robustness of the method. The error bars on the time 
delays were calculated using Monte Carlo simulations, i.e., re- 
distributing the magnitudes of the data points according to their 



photometric error bars. The width of the resulting time delay 
distributions gives us the Icr error bars. 

Because of microlensing we do not have access to the in- 
trinsic variations of the quasar. We represent microlensing in 
three of the light curves as a relative variation with respect to 
the fourth light curve, taken as a reference. We tested each of the 
four light curves in turn as a reference and kept the one that min- 
imised the residual microlensing variations to be modelled in the 
3 others. This is best verified with component B as a reference. 

With B as a reference light curve, we note that microlensing 
in C and D remains smooth and can therefore be modelled with a 
low-order polynomial drawn over the full length of the monitor- 
ing. However, A contains higher frequency variations that need 
to be accounted for in each season individually, as illustrated 
in Fig. [5] In doing this, we obtain fairly good fits to the light 
curves, as shown in the residual signal. To quantify the quality of 
these residuals, we applied the so-called one -sample runs test of 
randomness, a statistical test to estimate whether successive rea- 
lizations of a random variable are independent or not. In practice 
the test was applied to a sequence of residuals to decide whether 
a model is a good representation of the data. For most seasons 
in our curves the number of runs was between Icr and 3<x lower 
than the value expected for independent random residuals. Thus, 
although our microlensing model is not fully representative of 
the real signal, the deviations from the data points remain small. 

We tested the robustness of our curve- shifting method in sev- 
eral ways. First, we modelled the microlensing variations using 
polynomial fits of different orders. Second, we fitted these poly- 
nomials either across each individual season or across groups of 
seasons. Finally, we masked the seasons with the worst residual 
signal (Fig. [5]). All these changes had only a negligible impact 
on the time delay measurements. We note that this is not the 
case when considering only two or three seasons of data, which 
shows the importance of a long-term monitoring with good tem- 
poral sampling. 

4.3. Final results 

Our results are summarised i n Table |4| and are com pared with 
the previous measurements of Kochan ek et al. (2006 ), who have 
used pure a polynomial fit to the light curves and two sea- 
sons of monitoring. Using the same data but with our modified 
dispersion technique, we obtained very similar time delays as 
Kochane k et al.| ([2006 ), but larger error bars. We prefer keeping 
a minimum possible number of degrees of freedom (e.g., in the 
polynomial order used to represent microlensing), in accordance 
with the Occam's razor principle, even to the cost of apparently 
larger formal error bars. 

We also note that Kochane k et al. ( 2006 ) give their time de- 
lays with respect to A, which, with seven seasons of data, turns 
out to be the most affected by microlensing. As a consequence 
the error bars on these time delays are dominated by residual mi- 
crolensing rather than by statistical errors. The time delays used 
in the rest of our analysis are therefore measured relative to B. 

Finally, we used the mean values of our microlensing cor- 
rections to estimate the macrolensing R-band flux ratios between 
the four quasar images, assuming that no long-term microlensing 
affects the data. We found m B - m A = Am B A = 0-62 ± 0.04, 
Am BC = 0.05 ± 0.01, and Am BD = -0.16 ± 0.01, which is well 
compatible with the ratios measured at seven wavelengths by 
Mosquera et al. (2010). However, these authors report signifi- 
cant wavelength dependence of the image flux ratios, which led 
us not to use flux ratios as a constraint in the lens models. 
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Fig. 5. Lightcurves obtained with all four telescopes and shifted by time delays of A^ B a = 8.4 days, A^ B c = 7.8 days and A^ B d = -6.5 
days. The relative microlensing representations applied on curves A, C and D are shown as continuous curves with respect to the 
dashed blue line (see text). A fifth-order polynomial was used over the seven seasons to model microlensing on the quasar images 
C and D, while seven independent third-order polynomials were used for image A. The lower panels show residuals obtained by 
subtracting a simultaneous spline fit (grey) from the light curves. 



Table 4. Time delays for HE 0435-1223, with the same arrival order convention as Kochanek et al. (2006), i.e., D arrives last. 



Data 


Method 


A^AB 


A^AC 


A^AD A?BC 


A^bd 


A^cd 


SMARTS (seasons 1 and 2) 


Kochanek (2006) 


-8.0 + 0.8 


-2.1 + 0.8 


-14.4 + 0.8 






SMARTS (seasons 1 and 2) 


dispersion 


-8.8 + 2.4 


-2.0 + 2.7 


-14.7 + 2.0 6.8 + 2.7 


-5.9 + 1.7 


-12.7 + 2.5 


COSMOGRAIL (all seasons) 


dispersion 


-8.4 + 2.1 


-0.6 + 2.3 


-14.9 + 2.1 7.8 + 0.8 


-6.5 + 0.7 


-14.3 + 0.8 



5. Constraining the mass profile of the lensing 
galaxy 

The goal of the present section is to constrain the radial mass 
profile of the lensing galaxy as much as possible, which is the 
main source of uncertainty on the determination of Hq with the 
time delay method. One way of doing this is to use the infor- 
mation contained in the Einstein ring that is formed by the host 
galaxy of the lensed quasar ( [Suyu et al.|2010 2009 1 (Warren & 
|Dye||2003| ). This works well when a prominent Einstein ring is 
visible. Unfortunately, given the depth of the current HST im- 



ages of HE 0435-1223, the radial extent of the ring is too small 
to efficiently apply this technique. We propose instead to use in- 
formation on the dynamics and on the stellar mass of the lens, 
using deep optical spectroscopy. 



5.1. Stellar population and velocity dispersion of the lensing 
galaxy 

A deep VLT spectrum of the lensing galaxy is available from 
Eigenbrod et al. ( 2006). While the spectrum was originally used 
to measure the redshift of the galaxy, it turns out to be deep 
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Table 5. Model parameters of the lens potential well. For the variables used as parameters the ranges of values used are given 
between brackets. 



R eff 


8.44 kpc 


radius containing 50% of the observed light 


Re 


6.66 kpc 


Fmstem radius 


fb 


parameter 


baryonic fraction in the Einstein radius [0.05 - 0.5] 


r* 


5.3 kpc 


stellar component scaling radius 


M* 


M h f b /(l-f b ) 


stellar total mass 


7%r,max 


20 r* 


stellar truncation radius 


7dm 


parameter 


dark matter inner slope [0 - 2] 


r s 


parameter 


dark matter scaling radius [1, 2, 4, 8] x R E 


M h 


parameter 


dark matter total mass [4.8 x 10 11 - 9.1 x 10 12 ]M o 


^s,max 


10 r s 


dark matter truncation radius 



enough (< S/N >~ 20) to measure the stellar velocity disper- 
sion and the mass-to-light ratio. 

We analysed the data using full spectrum fitting with the 
ULySS package (Koleva et al. 2009). The method consists in 
fitting spectra against a grid of stellar population models con- 
volved by a line-of-sight velocity distribution. A single minimi- 
sation allowed us to determine the population parameters (age 
and metallicity) and the kinematics (redshift and velocity dis- 
persion). The stellar mass-to-light ratio was derived from the 
age and metallicity for the considered model. The simultaneous 
fit of the kinematical and stellar population parameters reduced 
the de generacies bet ween age, metallicity and velocity disper- 
sion QKoleva et al.|2008) . 

We performed the spectral fi t using the PEGASE-HR single 
stellar population models (SSP, Le Borgne et al. 2004 ), where 
the observed flux, F A , is modelled as follows: 



Fa = P n (A)x[L(v sys ,cr)( 

+ GmW). 



)S(t, [Fe/H],A)] 



(D 



The models were built using the Elodie.3.1 (Prug niel & 
Soubiran 2001; Prugn ieTet al.||2007] ) spectral library and the 
Kroupa pCroupa| |2001| ) and Salpeter ( |Salpeter| [T955] ) initial 



mass functions (IMF). L(y sys ,cr) is a Gaussian function of the 
systematic velocity, v^, and of the velocity dispersion, <x. 
S (t, [Fe/H], X) is the model for the SSP and depends on age and 
metallicity. P n is a polynomial of degree n, which models any 
residual uncertainty in the flux calibration and extinction cor- 
rection. In addition, the quasar spectrum might not be perfectly 
subtracted from the galaxy. To mimic this effect, we included 
an additive polynomial, Q m . The order of the additive and mul- 
tiplicative polynomials is the minimum required to provide an 
acceptable^ 2 , i.e., in our case (n = 10, p = 1). 

We obtained SSP-equivalent ages and metallicities of 
t ~ 3Gyr and [Fe/H] ~ 0.0 dex respectively. The cor- 
responding rest-frame B-band stellar mass-to-light ratio is 
M+/L B = 3.2+^ Mq/Lq^b using a Kroupa IMF and M+/L B = 
4.6+Q 7 Mq/Lq^b using a Salpeter IMF. The uncertainties in the 
age and metallicity were estimated via Monte Carlo simulations 
and propagated in the error in M*/L 5 . 

To compute the physical velocity dispersion we subtracted 
quadratically the instrumental broadening from the measured 
profile, neglecting the dispersion of the models since they are 
based on high-resolution templates. The instrumental broaden- 
ing was measured both from the PSF stars used to carry out the 
spatial deblending of the spectrum (Eigen brod et al.||2006| ) and 
from the lamp spectra. We obtained the rest-frame physical stel- 



lar velocity dispersion of the lensing galaxy: <x ap = 222 ± 34 km 
s" 1 in an aperture of 1", i.e., 5.7 kpc. 

5.2. Numerical integration of the Jeans equations 

In this section we model the lensing galaxy using a 3D spheri- 
cal potential well formed of two components, one for the stel- 
lar part of the mass and one for the dark matter halo. We then 
perform a numerical integration of the Jeans equations in 3D 
to predict a theoretical velocity dispersion and a total mass for 
the model. The assumption of spherical symmetry is sufficient 
for our purpose, as ill ustrated by the study of the lensed quasar 
MG 2016+112, where [Koopmans & Treu| ( |2002| ) introduced an 
anisotropy parameter and showed that it has almost no influence 
on the inferred mass slope. 

The luminous component of the model is a Hernquist profile 
( Hernquist|1990l ): 



P*(0) 



(r/r*)(l + r/r*) J 



(2) 



where p*(0) is the central density and r+ is a scale radius chosen 
so that the integrated mass in a cylinder of radius R q q (effective 
radius) is equal to half the total stellar mass M*. The profile has 
a maximum radius of r+ max = 20 r+. 

The dark matter halo is modell ed as a generali sed Navarro, 
Frenk & White (NFW) profile ( [Navarro et al.|1996| ): 



Ph(r) = 



Ph(0) 



(r/r s y DM (l + (r/rs) 2 ) (3_rz)M)/2 ' 



(3) 



where y DM is the inner slope of the profile, r s is the scaling radius 
and ph(0) is the central mass density. For y DM = 1 the model 
closely follows the standard NFW profile. Its total mass, Mh, is 
given in the truncation radius max = 10 r,. 

Following the usual convention, the integrated stellar and 
dark matter masses are related by the baryonic fraction fb, 



f b = 



M h + M+ 



(4) 



Because all integrations in this work were carried out nu- 
merically, fb can easily be computed in any aperture. We chose 
to compute it in the Einstein radius, which is also where lensing 
gives the most accurate mass measurement. The velocity dis- 
persion of the stellar component was computed by solving the 
second moment of the Jean s equation in spherical coordinates 
(Binney & Merrifield 1998 ). The velocity dispersion is then 



t(r) 



1 



P*(r) 



£ 



dr' Pir {r')d f <$>(r'), 



(5) 
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Fig. 6. Radial slope, Jdm, of the dark matter halo as a function 
of the baryonic fraction in the Einstein radius, R E . The colour 
code represents the value of^j eans for all acceptable models, i.e., 

^Jeans < ^ ( see text; ec l- ^he observational limits drawn on fb 
using our spectrum of tne lensing galaxy are indicated as shaded 
vertical regions. The blue region corresponds to a Kroupa IMF 
and the red regions corresponds to a Salpeter IMF. In this figure, 
all models have r s = 4 x R E = 25 kpc. 



where O(r) is the total gravitational potential. Equation ^ 
is solved numerically as follows. First, the density of each mass 
component (stars and dark halo) is sampled by N particles, using 
a Monte-Carlo meth od. Second, the poten tial is computed using 
a treecode method (Barnes & Hut 1986). For self-consistency, 
the density is computed by binning the particles in spherical 
shells. Then, a velocity is allocated to each particle at a distance r 
from the galaxy center, following a Gaussian distribution of vari- 
ance <r*(r). Finally, the velocity dispersion of the model, <xj ean s> 
is computed numerically, by integrating all the particle velocities 
in an aperture that matches exactly the slit used to carry out the 
observations. 

The advantages of the numerical representation of the lens 
models are multiple. They allow us (i) to compute velocity 
dispersions for any combination of density profiles, even non- 
parametric ones, (ii) to account for the truncation radius of the 
halo, (iii) to test the dynamical equilibrium of the system by 
evolving the model with time. Computing the velocity disper- 
sion as well as the total mass in a cylinder along the line of sight 
of the observer is then straightforward. 

In all calculations we fixed the effective radius to the one ob- 
served for the lensing galaxy, i.e., R q q - 8.44 kpc. The remaining 
parameter space to explore is composed of the halo scale radius, 
r s , the slope of the profile, Jdm, and the baryonic fraction, fy, 
within the Einstein radius. 



5.3. Applying the dynamical and stellar population 
constraints 

Using our HST photometry, we measured the total rest-frame 
B-band galaxy luminosity by converting its total H-band mag- 
nitude, rrifrnw = 16.20. Using a k-correction of 1.148, and a 




0.4 0.6 

Baryonic Fraction in Einstein Radius 

Fig. 7. Effect of change in r s , the scale radius of the dark matter 
halo, on the y DM - fb relation. The contours show the region 
containing 68% of the models. The blue area has r s = 25 kpc 
(our adopted choice for the final modelling), the green region 
has r s = 15 kpc, and the red region is for models with small 
scale radii, r s = 5 kpc. 



galactic extinction^ E(B - V) = 0.059 ( [Schlegel et al.|l998| , 
we find L B = 1.04 xlO 11 L fi?0 . For a given galaxy model we can 
therefore compute the total stellar mass-to-light ratio, M+/L B , as 
well as the baryonic fraction that we can compare with the ob- 
served ones. This requires a choice of IMF. For a Salpeter IMF 
we measured ft = 0.65^ }^, while for a Kroupa IMF we mea- 
sured fb = 0.45^qq^. In our computation of the baryonic fraction, 
the total mass is the mass in the Einstein radius, M(< R E ). 

We compared the model properties to the data by computing 
a chi-square,^j eans . This^j eans includes the measured lens veloc- 
ity dispersion and the total mass in the Einstein radius, M(< R E ) 
= (3.16+0.31) xlO 11 M . We allowed for a 10% error on the 
total mass in the Einstein radius to account for the weak depen- 
dence of this mass upon the choice of a lens model, which led 
to 



2 

X Jeans 



/ cr Jeans -222 \ 2 + ( M(< R E ) - 3.16 x 10 1 



34 



0.31 



(6) 



We used our numerical integration of the Jeans equations to 
sample the y DM - fb plane, for a given choice of the scale radius, 
r s , of the dark matter halo. Fig. shows the value of X^ Qms f° r 
960 000 galaxy models with r s = 4 x R E = 25 kpc. For a given 
choice of IMF a precise measurement of the baryonic fraction, 
fb, can constrain the dark matter slope jdm well, at least in the 
case of a Kroupa IMF. 

Fig. [6] is drawn for r s = 4 x R E = 25 kpc, which gives the 
best ^j eans values. However, models with Afj eans < 1 are also 
obtained with other values of r s . We therefore test the influ- 
ence of a change in r s by building models for a range of scale 
radii, r s = 5, 15,30 kpc (i.e., r s = 0.75x, 2.25x, 4.50 x R E ). 
Fig. [7] summarises the results and shows no influence of r s on 
the y DM - fb relation unless r s becomes unrealistically small. 

The Hernquist profile in our model also has a scale parame- 
ter, r*, which is fixed to the observed value of the effective ra- 



1 A galactic extinction calculator is available at this address: 
http ://ned w w w. ipac . caltech . edu/f orms/c alculator. html 
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Fig. 9. Distribution of our lens models as a function of the dark matter slope and Hq. The colour code gives the value of Xj ot (see 
text; eq.|7j). Only the models with^ ot < 1 are shown. The results in the left panel are for a Salpeter IMF, and for a Kroupa IMF in 
the right panel. The lensing;^ 2 itself is not included in the figure. However, it is systematically lower for galaxies with Kroupa IMFs 
than for galaxies with Salpeter IMFs. 
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Fig. 10. Same as Fig. [9] but the colour code now gives the slope of the total (dark+luminous) mass profile, y' . A logarithmic slope 
of Y =2 corresponds to an isothermal profile, shown in green. 



dius, 7? e ff = 8.44 kpc. Although we have HST images that allow 
us to measure and although our error bars on this parameter 
are small, we note that studies of lensing galaxies by different au- 
thors often find different 7? e ff values from the same data. This is 
because different methods are used to subtract the quasar images 
and to mask (or not) the Einstein ring. We therefore computed 
models for a range of different and show the influence on 
the y DM - f b plane in Fig. [8} for # eff = 4.22, 8.44, 12.88 kpc. 
This reflects an (unlikely high) 50% change in R q q with respect 
to the actual measurement and translates into a 10-20% shift of 
the y DM - fb relation. This is comparable to the width of the 
coloured areas in Fig. [8] Because our true measurement error on 
7? e ff is certainly much better than 50%, we conclude that this pa- 
rameter does not affect the models much, given the accuracy of 
the observations. 



6. Towards H with HE 0435-1223 



The gravlens software ( |Keeton|[2004| ) was used in combina- 
tion with the dynamical and stellar population constraints to re- 
produce the lensing configuration of the quasar images and to 
attempt converting the time delays into Hq. We adopted a total 
potential well composed of the main lensing galaxy, the nearby 
galaxy G22 (see Morga n et al.|20 05) plus an external shear. 

The main galaxy is composed of a projected Hernquist + 
cuspy halo model identical to those in eqs. 1 & 2. 

For the Hernquist profile, we fixed the ellipticity and PA of 
the lensing galaxy to the observed ones. 

We then estimated how well each of the models described in 
Section [5^2| reproduces the observed image configuration. In do- 
ing this, we allowed only the external shear (y, # r ), the Einstein 
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Fig. 8. Effet of a change in 7? ef f on the y DM - fb relation. The 
contours show the region containing 68% of the models. The 
value adopted in this paper is 7? e ff = 8.44 kpc, as measured from 
HST images. This corresponds to the area in red. The models in 
the blue region have 7? e ff = 4.22 kpc, and the region in green has 
# eff = 12.88 kpc. 



radius of G22 and Ho to vary. The potential well of G22 was as- 
sumed to lie at its observed position and was modelled as a sin- 
gular isothermal sphere (SIS). We assumed a conservative valu e 
ofR E (G22) < a/4, following the results of |Morgan et aLlpOQ5) , 
who have found R E (G22) =0718. 

We show in Fig. [9] the value of Hq for each lens model as 
a function of its dark matter slope, Jdm- The colour code in the 



figure corresponds to the value ofx\ oV where the baryonic frac- 



tion, f b ± cr(f b ), 
calculation: 



in the Einstein radius is now included in the 



%Tot % Jeans ~'~ 



f b (model) - /j,(obs) 

oifb) 



(J) 



Including fb is justified by Fig.[6| showing that different va- 
lues of fb select different lens models. We display our results for 
the two most common IMFs in use in stellar population mod- 
elling, the Salpeter and the Kroupa IMFs. 

The points define a x\ oi surface with a clear valley which 
minimum indicates the best dark matter slopes for each IMF. 
These are yz)M(Sal) ~ 1.15 and ^^(Kro) ~ 1.54 for the Salpeter 
and the Kroupa IMFs, respectively. Each model shown in Fig. [9] 
is also required to display a good lensing chi- square,;^, after fit- 
ting of the quasar image positions with gravlens. The values of 
x\ are systematically lower for the lensing galaxies with Kroupa 
IMFs. The present lensing and dynamical work therefore favours 
a lensing galaxy with a Kroupa IMF, as also found by Cappellari 



et al. (2006 ) fr om 3D spectroscopy of early-type galaxies and by 



(2008) using the SLACS sample of strong lenses 



Ferreras et al.[ 
( [Bolton et al.|2008| ). 

The points in Fig.[9]all have^j eans < 1, which prevents us, for 
now, from giving a value for Ho given the observational uncer- 
tainty on fb and <x ap . The measurement errors on both parameters 
will, however, easily improve with deeper and higher resolution 
spectroscopy of the lens. 

With the current observational constraints we rely on pre- 
vio us work done on the tot al mass slope, y\ of lensing galax- 
ies. Koopmans et al. ( 2009| ) have measured the probability dis- 



tribution function of y' in the SLACS sample of strong lenses. 
They have found (y'> = 2.085^^, with an intrinsic sp read of 
o~(y r ) = 0.2 0, also confirmed in a more recent study bv Auger 



et al. 



2010). The models with the best^ ot values in 



ly by_i 
Fig-g 



corre- 



spond to a total slope of y' - 2. 1 ± 0. 1 independent of the choice 
of an IMF. This is well within the rms limits of Koopmans et al. 
(2009) and leads to 57 < H < 71 km s" 1 Mpc A with all models 
equiprobable within this range. 

In our analysis we modelled the environment of the lens- 
ing galaxy as a SIS that represents galaxy G22, plus an exter- 
nal shear with a PA of the order of -15°, and amplitudes in the 
range 0.06 < y < 0.08. However, we did not explicitly account 
for the fact that H E 0435-1223 lies within a group of galaxies 
(Momcheva 2009). The unkown convergence, /c, associated with 
the group leads to an overestimate of Hq by a factor (1 - k)~ 1 , 
meaning that the range of acceptable Hq values would decrease 
further if we have underestimated the convergence caused by the 
group. Current imaging data suggest, however, that the line of 
sight is i n fact sightly under-d ense compared with other lensed 
systems (Fassnacht et al. 2011). However, the missing conver- 
gence is about k ~ 0.01 - 0.02, leading to at most a 2% change 
(upwards) in Hq. The above statement should, however, be con- 
sidered with care because the effect of a group is poorly appro- 
ximated by a simple convergence term. Explicit modelling of the 
group halo is likely needed to properly account for the modifi- 
cations induced on the main lens potential. Deep X-ray and/or 
optical integral field spectroscopy may turn out to be very useful 
as well in determining the centroid and mass of the group that 
contains the galaxy lensing HE 0435-1223. 



7. Conclusion 

We presented seven years of optical monitoring for the four 
lensed quasar images of HE 0435-1223. We found that the time 
delays are better expressed with respect to component B, which 
is the least affected by stellar microlensing in the lensing galaxy. 
The formal error bars on the time delays are between 5% and 
10% depending on the component, which is remarkable given 
the very short time delays involved. In addition, the delays are 
robust against different tests performed on the data, including 
removal of subsets of data and monte carlo simulations. These 
tests are possible only with very long light curves, as provided 
by COSMOGRAIL. Most past lens monitorings have two to 
three seasons and much coarser temporal sampling. The delays 
are also independent of the way the microlensing variations are 
modelled. Given the short length of the delays, additionally im- 
proving the error bars will only be possible by increasing the 
temporal sampling of the curves, i.e., by merging all existing 
data on HE 0435-1223 taken by different groups over the years. 

We introduced a method to convert the time delays into Ho 
purely based on external constraints on the radial mass slope of 
the lens. These constraints come from deep optical spectroscopy 
of the lens, which allowed us to measure its velocity dispersion 
and its stellar mass-to-light ratio. The present paper describes 
our approach, but the current observations of the lensing galaxy 
so far lead to a broad range of values for Hq. 



Our methodology complement that presented by Suyu et al. 
( 2009, 2010) well. Combining the two in the future will require 
follow-up observations such as (i) deep HST imaging to map the 
Einstein ring with high signal-to-noise, (ii) deep high resolution 
spectroscopy of the lens over a broad spectral range to narrow 
down the uncertainty on its velocity dispersion and stellar mass, 
(iii) integral field spectroscopy to measure all redshifts within 
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30-60" around the lens, and (iv) X-ray imaging to pinpoint mas- 
sive groups along the line of sight. 

If all these observations are taken for HE 0435-1223 as well 
as for a few lenses with well measured time delays, the cost in 
terms of follow-up observations will still remain very modest 
compared with other existing methods to measure Hq, such as 
Cepheids and supernovae. 
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